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Abstract 

Background: Disease is a major factor driving the evolution of many organisms. In honey bees, selection for social 
behavioral responses is the primary adaptive process facilitating disease resistance. One such process, hygienic 
behavior, enables bees to resist multiple diseases, including the damaging parasitic mite Varroa destructor. The 
genetic elements and biochemical factors that drive the expression of these adaptations are currently unknown. 
Proteomics provides a tool to identify proteins that control behavioral processes, and these proteins can be used 
as biomarkers to aid identification of disease tolerant colonies. 

Results: We sampled a large cohort of commercial queen lineages, recording overall mite infestation, hygiene, and 
the specific hygienic response to V. destructor. We performed proteome-wide correlation analyses in larval 
integument and adult antennae, identifying several proteins highly predictive of behavior and reduced hive 
infestation. In the larva, response to wounding was identified as a key adaptive process leading to reduced 
infestation, and chitin biosynthesis and immune responses appear to represent important disease resistant 
adaptations. The speed of hygienic behavior may be underpinned by changes in the antenna proteome, and 
chemosensory and neurological processes could also provide specificity for detection of V. destructor in antennae. 

Conclusions: Our results provide, for the first time, some insight into how complex behavioural adaptations 
manifest in the proteome of honey bees. The most important biochemical correlations provide clues as to the 
underlying molecular mechanisms of social and innate immunity of honey bees. Such changes are indicative of 
potential divergence in processes controlling the hive-worker maturation. 
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Background 

Social insects such as the honey bee (Apis mellifera L.) 
derive great benefit from living in tight-knit groups that 
enable greater efficiencies in brood care, foraging and 
defense against predation. However, the high population 
densities and relatedness of individuals leave colonies sus- 
ceptible to emerging infectious diseases [1]. Varroa 
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destructor, an ectoparasitic mite of the honey bee [2] 
causes varroasis, which is a leading contributor to ongoing 
colony losses in commercial apiculture worldwide [3]. 
V. destructor feeds on the hemolymph of larval and adult 
bees, inflicting nutritional stress and immune suppression, 
as well as acting as a major vector for viral pathogen 
transmission [4]. 

In solitary insects, cellular or humoral-based defenses 
provide the only known system for immunity, but 
A. mellifera's genome reveals that while honey bees con- 
tain these systems for immunity, the number of immu- 
nity genes is lower than that of solitary insects such as 
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flies, moths and mosquitoes [5]. As an apparent com- 
pensation for this, social insects have evolved collective 
systems of behavior that provide defenses against disease 
and parasitism. Two related behaviors, hygienic behavior 
(HB) and Varroa sensitive hygiene (VSH), are highly 
variable among A. mellifera colonies and are seen as 
important traits in the development of disease and mite- 
resistant stock. HB is a well-documented protective 
behavior that involves nurse-aged worker bees uncap- 
ping brood cells and removing parasitized or diseased 
pupae [6]. VSH is less well-understood but it encom- 
passes a suite of behaviors that ultimately suppress mite 
reproduction by uncapping and/or removing mite- 
infested pupae from sealed brood resulting in a high 
proportion of non-reproductive mites in the brood that 
remains [7,8]. HB and VSH can be quantified using field 
assays and are heritable so, while both are now used in 
the selective breeding of Varroa-resistant bees [9,10], 
the genetic and biochemical mechanisms that drive 
them are poorly resolved. 

To date, most selective breeding in commercial apicul- 
ture focuses on traits such as honey production, color, 
gentleness, winter survival or other economic parameters. 
When combined with continual dilution of the gene pool 
through importation of susceptible stock, these selections 
limit host adaptation to pathogens. In order to improve 
disease and mite tolerance, field assays for HB and VSH 
must be incorporated into the stock selection process 
[11,12]; however, these assays are resource intensive, lack 
sensitivity and may require closed breeding [13], limiting 
their suitability for widespread application. To support the 
creation of novel assays, a molecular-level mechanistic 
understanding of resistance traits is seen as a promising 
avenue to support commercial breeding and disease pre- 
vention through marker-assisted selection (MAS) [14]. To 
date, low-resolution microsatellite-based quantitative trait 
loci (QTL) for HB have been reported [15], as have some 
of the biochemical consequences to the host of infection 
by V. destructor and associated viruses [16,17]. Transcrip- 
tome changes in A mellifera and in Varroa's natural host 
A. cerana also pinpoint subtle changes in transcript 
expression for components responsible for neuronal rewir- 
ing, olfaction, metabolism and aspects of social behavior 
that may be critical components driving mechanisms of 
Varroa tolerance [18,19]. 

All the molecular investigations of HB and VSH have 
used well-controlled colonies or individual bees without 
examining the natural variation and distribution of both 
the traits and their molecular components. Thus, here we 
tested the hypothesis that inter-colony variation in disease 
resistance parameters is reflected by changes in the expres- 
sion of specific proteins. Sampling from a large cohort of 
colonies, we measured the relative abundance of approxi- 
mately 1,200 proteins from two bee tissues involved in 



interactions with the pathogens and correlated these with 
estimates for active bee behavioral phenotypes for HB and 
VSH, as well as host-pathogen population dynamics. 
Through meta-analysis of these data with other available 
information, proteins and biochemical processes most 
likely to be responsible for the observed disease resistance 
traits were identified. 

Results 

Behavioral phenotypes and mite-bee dynamics 

To assess the expression of disease tolerant behaviors, 
colony-level measurements of various metrics of HB, 
VSH and mite-host population dynamics were made 
(Figure 1), including phoretic infestation (PH, mites on 
adult bees), natural mite drop (ND, estimate of mite 
death rate), and levels of brood infestation (BI). HB was 
estimated by observing the proportion of a defined num- 
ber of freeze-killed, sealed brood cells that bees first 
uncapped (uncapped, U) and then removed (removed, R) 
at 24 and 48 hours (Figure 2a). The hygienic response to 
freeze-killed brood was time dependent, with widely dis- 
tributed levels of HB at 24 hours and the majority of 
hives achieving the accepted 95% threshold for the pro- 
portion of removed cells at 48 hours [11]. Because colony 
scores in the other measured parameters were distributed 
roughly similarly to HB at 24 hours as indicated in the 
kernel density plot of Figure 2b, we asked which of them 
were independent measures and which were interrelated. 
Pearson's product-moment correlation (PPMC) coeffi- 
cient revealed a statistically significant negative pair -wise 
dependence between estimates for hygiene and mite 
infestation dynamics (PH, ND) (Figure 2b). Maximum 
dependence was observed between cell removal and PH 
(r = -0.54, P = 0.0007), confirming that colonies better 
able to detect and uncap cells with affected brood are 
able to reduce adult infestations more efficiently. Obser- 
ving differences in density for temporal aspects of 
hygiene and interactions between disease-tolerant beha- 
viors and infestation indicated we were detecting natural 
variation in the speed of HB removal that directly influ- 
enced the colonies' tolerance to Varroa mites. 

Correlation of protein expression and behavior 

While most of the parameters discussed above are known 
to have a genetic basis, all must ultimately manifest as the 
result of changes in protein expression and/or activity. 
To explore potential mechanisms underpinning natural 
variation in Varroa tolerance across these colonies, we 
examined the protein expression profiles of two tissues 
that play a critical role in the bee-Varroa interaction: 
antennae of brood-nest workers (that is, mostly nurse 
bees) and the integument from fifth instar worker larvae. 
Antennae were used because they are adult bees' primary 
sensory organs and many of the behaviors evaluated here 
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Figure 1 Diagram depicting honey bee disease tolerant traits and infestation dynamics (a) Hygienic behavior (HB) is composed of two 
component behaviors, 'uncapping' (uncapped, U) which involves the opening of the cell containing a dead pupa and 'removal' (removed, R) 
which involves the removal of the dead pupa from the cell after uncapping has occurred. These behaviors are not always performed by the 
same bee. HB was recorded over 24 hour (rapid) and 48 hour (slow) periods, (b) Varroa sensitive hygiene (VSH) was defined by determining the 
proportion of farroa-infested cells in which no reproductively viable Varroa mite daughters were produced. Increases in this measure infer that 
greater proportions of mites remaining in the brood have had their reproduction suppressed because of infertility, death, the production of only 
males, or have had their reproduction delayed preventing sexual maturation of females, (c) Brood infestation (Bl) is the percentage of brood cells 
infested by one or more mites regardless of the mite's reproductive status, (d) Phoretic infestation (PI) is an estimate of the density of mite 
phoresy on adult bees, and (e) natural drop (ND) is a normalized measure of the number of mites falling from the adult bees onto an adhesive 
board on the bottom board of colonies. 



involve bees being able to sense the presence of either the 
pathogen itself or a damaged/diseased nest-mate. Integu- 
ment was chosen because it is the initial physical barrier 
to Varroa when they feed on larvae and as such the 
innate processes found here may be critical components 
in the response of hygienic adults and provide direct 
innate mechanisms of tolerance. It is possible that 
changes in the composition of larval proteins or the 
metabolites produced by these proteins during infection 
may trigger HB or VSH in adult nurses. 

Using liquid chromatography-tandem mass spectrome- 
try (LC-MS/MS) to analyze three independent samples of 
each tissue, we constructed protein expression profiles 



for approximately 1,200 proteins across all colonies as 
described previously [20] . By centering and standardizing 
across labels and colonies, the relative expression ratios 
from individual LC-MS/MS experiments are converted 
into a roughly normalized distribution of protein effect, 
representing the expression level of each protein in each 
colony relative to the population average. These variables 
(response variable) were then regressed against the beha- 
vior and infestation estimates (predictor variable) mea- 
sured for that colony. The direction of each regression 
was determined by the sign of the estimated regression 
coefficient and the significance of that effect was accessed 
using a mixed linear model with probability cut-off at 
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Figure 2 Distribution of and correlation between behavioral traits, (a) Box and whisker, and density plot for disease tolerant traits for 
hygiene at different times after the introduction of freeze killed brood and the measure for VSH. (b) Correlation matrix showing the correlation 
score and significance for estimates of rapid hygienic behavior and VSH with measures of infestation. VSH, Varros sensitive hygiene. 
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of proteome correlation screen for antennal and larval tissues. Bar plot giving the number of proteins found correlating 
different cut-off's for statistical confidence. 



Q <.2 adjusting for multiple comparisons or later P < .05 
for explorative data analysis (see Figure 3). 

Several proteins are highly significant predictors of 
resistance to Varroa mite infestation 

To adjust significance levels to account for the multiple- 
testing hypothesis, proteins were filtered using Q <.2 
cut-off; for HB one antennal protein and five larval pro- 
teins survived this additional filter (Tables 1, 2). In the 
antennae, the hypothetical protein 'LOC552009' of 
unknown function correlated with HB at 48 hours (Q = 
0.09) for both 'uncapped' and 'removed' behaviors. 
Sequence analysis revealed that LOC552009 contains a 
conserved domain similar to the mammalian protein 
lipid transport protein Apolipoprotein O (ApoO) [21]. 
Figure 4a, b shows the added variable plot for this pro- 
tein correlating with HB (removed 48 hours), peptides 
identified and protein sequence containing the con- 
served domain for ApoE. 

In larvae, several more candidate proteins were identi- 
fied as strong positive and negative predictors for HB 



(Table 2), suggesting that events in the larvae may be 
able to influence HB of the adult. Further correlation 
analysis of mite infestation/fertility measures (PH, BI, ND 
and VSH) identified the hemocyte protein-glutamine 
gamma-glutamyltransferase-like (a putative transglutami- 
nase) as highly significant (Q = 0.02) and positively corre- 
lated with ND (Figure 4c). 

To increase the specificity of our measures for infestation 
dynamics, we next calculated the ratio of mites observed 
phoretically to those found in brood cells (PH/BI). This 
adjustment enabled quantitation of the relationship 
between two important stages in the mite life cycle, where 
long phoretic phases may be indicative of poor reproduc- 
tive success and the influence of adult bee behavior or lar- 
val attractiveness. After adjustment, several proteins (nine 
in antenna, four in larva) were highly significant (Q <.2) in 
both tissues. Importantly, the adjusted metric was also 
highly correlated with ND (Pearsons product correlation: t 
= 5.2211, df = 36, P = 7.633e-06) and in larva correlated 
with increased significance with the protein Tg (Q = 0.01), 
supporting the role of Tg as an accurate measure of Varroa 
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Table 1 Name, expression and function of larval proteins found correlating with traits measured for resistance to 
Varroa infestation. 



Protein Name 


R48 U48 PH/BI 


ND Function 


similar to CG5903-PA 


I 1 


CDD:ApoE (Apolipoprotein, lipid transport) 


PREDICTED: polyadenylate-binding protein 1 -like isoform 1 


1 


mRNA surveillance pathway 


PREDICTED: 60S ribosomal protein L14 


1 


Translation 


PREDICTED: alpha-tocopherol transfer protein-like 


I 


Vitamin E distribution 


PREDICTED: alcohol dehydrogenase [NADP+] A-like 


1 


Alcohol metabolism 


phenoloxidase subunit A3 


t 


Monooxygenase, phenoloxidase activity 


PREDICTED: dehydrogenase/reductase SDR family member 11-like 


t 


Oxidoreductase activity (Metabolism) 


PREDICTED: v-type proton ATPase subunit F Mike 


T 


Oxidative phosphorylation (Electrom transport) 


PREDICTED: calcyphosin-like protein-like 


r 


Intracellular signal transduction 


PREDICTED: sorbitol dehydrogenase-like isoform 2 


r 


Fructose and mannose metabolism 



Bl, brood infestation; ND, natural drop, PH phoretic infestation; R, removed; U, uncapped. 



resistance (Figure 4, Tables 1, 2). Tg also provides the 
clearest link to phenotype. In insects, Tg is the primary 
protein facilitating the crosslinking of the clotting factors, 
hemolectin and Eig71Ee. Tg drives clot formation, an 
important defensive process against ectoparasites such as 
the Varroa mite. Of the proteins correlating in the antenna, 
several form part of diverse metabolic pathways and are 
impossible to link to function at this point. However, calcy- 
phosin-like protein and phenoloxidase subunit A3 are 
worthy of discussion. Calcyphosin (Q = 0.02) has been 
found expressed in olfactory cells in lobster and with a 
putative role in signal transduction in sensory cells. Pheno- 
loxidases (Q = 0.1) are a critical part of insect immunity 
during the pathogen encapsulation response, but are also 
important in the ecdysteroid-dependent processes linked 
to polyphenism and caste differentiation. 

Neuronal proteins underpin hygienic behavior and VSH 
in antennae 

The proteins discussed above that survive correction for 
multiple hypothesis testing should be excellent predictors 
of HB, perhaps even usable in marker-assisted selective 
breeding. To find so many highly significant proteins in a 
completely out-bred population is remarkable but the 



requirements that they must pass mean that it may be 
too restrictive a dataset to understand fully some of the 
molecular mechanisms underlying the relevant behaviors. 
To this end, we expanded the analysis to discover pro- 
cesses with mechanistic relevance for HB and VSH. 
Those proteins with a significant correlation (P < .05) to 
one or more behavioral traits were explored using onto- 
logical classifications provided by the honey bee refseq 
entry in the National Center for Biotechnology Informa- 
tion (NCBI) or the flybase homolog. The most significant 
of these enrichments (P < .05, Table 3) largely tracked 
with the altered distributions for estimates of HB, sug- 
gesting some of the molecular mechanisms that may reg- 
ulate this behavior. The set of proteins highly correlated 
with rapid HB (>95% removed by 24 hours) was particu- 
larly enriched for proteins involved in 'sensory develop- 
ment'. At 24 hours, both uncapping and removal traits 
correlated with the up-regulation of the secretory pro- 
teins windbeutel, amphiphysin (Amph) and [RefSeq: 
CG6259] which encodes the homolog to human CHMP5 
protein. Proteins down regulated in rapid hygienic bees 
were ankyrin 2, laminin A, Zasp (Z band alternatively 
spliced PDZ-motif protein) and fasciclin 1 (Fasl) all 
involved in 'cell adhesion'. Proteins correlating with both 



Table 2 Name, expression and function of antennal proteins found correlating with traits measured for resistance to 
Varroa infestation. 



Protein name 



R48 U48 PH/BI ND Function 



eukaryotic translation initiation factor 3 subunit A II" 

rab proteins geranylgeranyltransferase component A 1 II" 
argininosuccinate synthase-like | 

glucose-6-phosphate 1-epimerase-like f f 

UDP-glucose:glycoprotein glucosyltransferase f 
mitochondrial import inner membrane translocase subunit Tim8 - - J, 

60S ribosomal protein L10 isoform 1 - - | 

26S proteasome non-ATPase regulatory subunit 2 " " T 
hemocyte protein-glutamine gamma-glutamyltransferase f 



Mitotic spindle elongation; mitotic spindle organization 
Sensory transduction, intracellular transport 
Alanine, aspartate and glutamate metabolism 
Glycolysis/Gluconeogenesis 
Protein glycosylation 

Protein import into mitochondrial inner membrane 

n/a 

Regulatory subunit of the 26 proteasome 
Peptide cross-linking, emolymph coagulation 



Bl, brood infestation; ND, natural drop, PH phoretic infestation; R, removed; U, uncapped. 
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(a) 



>gi | 32 8780825 | ref I XP_62 4 3 92 . 3 _PREDICTED : _hypothetical_protein 
_LOC552 00 9_isof orm_l_ [ Apis_mellif era ] 

mqrtlffkkflmpcglcaavpamkppipeehpapcsneiqgkklikpselpiysiddgyt 
kqmpciqypsivednirkirqtvseikltidkisrdisssleslkfi sdylqdqanlmpr 
igavgvgglsglilslrggifkrlmytttgaaivgcvcfpketketvntmehygnvsynf 
iygvkpgdnkkeisfnefplvksvleseyfrmlvqffeqktndtptttdwltdttakte 
inkkk! 



ApoO conserved donmain 



Peptides identified 



(c) 



>gi I 66510575 I ref I XP_392939. 2 I PREDICTED: hemocyte 
protein-glutamine gamma-glutamyltransf erase-like [Apis 
raellifera] 

MSREEPLWEMVYLYEKENAKLHHTINYELVHLDPPAAVLRRGQSFHIALRFNREYIDE 
IDIVRLLFSFGPNPNVLRGTRGVNTITNRDSYLTDLEAWGVRLIGVSGVDLSAEVRSPV 
DSPVGMWQLNIETTIVGSKRSPNTYNYDKDIYLLFNPWLKEDLVYMEDEQLLDEYILND 
VGKIWVGAWGSARGREWIFGQFDAYVLKACQLLLERSGIKANSRGDPIQMCRAISRIVN 
SNDDKGWTGRWDGDYQDGTAPAAWTGSVPILEQFLETGESVKYGQCWVFAGWTTVCR 
ALGIPSRWSNLVSAHDANASLSVDRYYSKENEELEYDPNNVEGEDSIWNYHVWNDVWM 
ARPDLPKGYGGWQAIDSTPQEPSEGVYQCGPASVEAIKQGWGYNYDVTFMLASVNADL 
MRWIEDPDSEMGFRKIDCNKYHIGRMILTKAPWVYDPNGDRDREDITSLYKAKEGTELE 
RLTLYRAVRSTELAKRFYSLPSPAKEDVEFDLVDIERVNIGEPFAVIVNIKNKSNEKRT 
IQAILSAGSVYYKGIKAYLVKRASGDFVLEPYASEQLRLTITVDDYLDKLVEYCNMKLY 
SIATWETKQTWADEDDFQVLKPNIVVKIDGEPTVGKPSIISLRFKNPLQRVLTDCKFN 
YAGPGLTRNKTLAFRDVDPEEDVYVEHQLIPQKAGSQKIIATFTSKELVDVTGSAVIDV 
LDMDE 



Domains 



Transglut (pfam00868) 
Transglut_core (cl10469) 
Transglut_C ( pfam00927) 



Peptides identified 




-20 -10 o 
Removed 48 h 




0.000 0.001 

Natural drop 



Figure 4 Protein sequence, peptide identification coverage, and conserved domains for two proteins with highly significant 
correlation with Varroa tolerant traits, (a-b) Sequence encodes a probable ApoO type protein found to correlate with hygienic behavior at 
48 hours, (c-d) A transglutaminase (Tg) correlates with natural drops, an estimate of mite death rate in the colony. ApoO, apolipoprotein O. 



elements of HB at 48 hours (>95% open/cleaned 48 h 
ours) were enriched in the ontology 'mitochondrial inner 
membrane'. This ontology was also identified for the 
rapid uncapping behavior and corresponds to reduction 
in primary metabolic pathways. 

VSH and HB correlated strongly with some of the same 
proteins, even though there was no apparent 



interdependence between VSH and HB (Figure 5). Pro- 
teins that correlated with VSH and HB included Fasl, 
which was negatively correlated with rapid HB and VSH, 
while Amph and helicase 25E (Hel25E) levels were signif- 
icant positive predictors of both VSH and rapid HB 
(Figure 5). In addition, several unique proteins with roles 
in synaptic function correlated only with VSH. VSH was 



Table 3 Summary of gene enrichments for larval proteins correlating with anti-parasitic traits. 

Behavior Expression Enriched ontologies 



R24 



U24 



R48 



U48 



VSH 



V 
413 

1 1 1 

414 



fl2 
47 
f14 
120 



Macromolecular localization (wbl, Amph, Hel25E, CG6259), protein localization, sensory organ development (wbl, Amph) 

tissue development, nerve growth factors (Ank2, LanA), axon guidance, cell adhesion (Ank2, LanA, Fas1, Zasp52, RpS25) 

Macromolecular localization, protein transport (wbl, Amph, Hel25E, CG6259), microtubule cytoskeletal organisation, cell cycle 
(RpS30, RpL12, Hel25E, Rab11) 

Cellular respiration (Gdh, ND42 CG5703 Cg6463), mitochodrial electron transport (Gdh, ND42 CG5703) 
Translation, mitotic spindle organization (Aats-cys, RpL21, Cctgamma) 
Mitochondrial inner membrane (sesB, ND42, CG5703, 1(2)06225) 
No enrichment 

Mitochondrial inner membrane (ND42, CG5703) 
Ubiquitin-dependent protein catabolic process (Rad23, Prosalpha7) 
Cytokinesis (Gammacop, SNAP) 



R, removed; U, uncapped; VSH, Varroa sensitive hygiene. 
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wbl 
Hel25E 
CHMP5 
Amph 
gammaCop 
Snap 
mub 
RpS25 
Fas1 
Zasp52 
Ank2 
LanA 




Function 



Macromolecular 
localisation 



Cytokinesis 



Cell adhesion 



T- Value 



1 



Figure 5 Heat map of t values for proteins with involvement in 
neuronal/synaptic function in antennae found to correlate 
significantly with hygienic behavior and VSH. VSH, Varroa 
sensitive hygiene. 



negatively correlated with expression of the mushroom 
body protein (Mub), soluble NSF attachment protein 
(SNAP) and gammaCOP. Mub is involved in temperature 
preference in Drosophila, [22] while SNAP is a key pre- 
synaptic protein mediating synaptic vesicle fusion and 
gammaCOP is involved in vesicle trafficking at synapses 
and other vesicle sorting pathways [23] . The protein with 
greatest change in expression with respect to VSH mea- 
surements was [RefSeq:LOC412768], a poorly annotated 
member of the take-out/juvenile hormone binding pro- 
tein (To/JHBP) superfamily involved in chemoreception. 

Larval proteins with cuticular and immune function 
correlate with HB and VSH 

HB and VSH are considered to be specific behavioral 
adaptations of the adult honey bee to diseased brood. 
However, as our molecular understanding of these traits 
improves, it is possible that factors expressed within the 
larva in response to the pathogens may influence the 
manifestation of the relevant behaviors. Protein/behavior 
correlations of larval proteins revealed enrichment of 



chitin-based cuticle structural proteins, particularly cuti- 
cular proteins 3 and 13 (CRP3 and CRP13; Table 4). 
Larval expression of the peptidyl-amino acid modifying 
enzymes Cafl and [RefSeq:CG6370] was also positively 
correlated with HB. CG6370 encodes a dolichyl-dipho- 
sphooligosaccharide-protein glycosyltransferase involved 
in N-glycan biosynthesis in the lumen of the endoplasmic 
reticulum and, interestingly, four other lumenal proteins 
were also positively correlated. The larval protein/beha- 
vior correlation with the steepest slope was [RefSeq: 
CG1318-PA], a beta-N-acetyl-D-hexosaminidase with 
broad substrate specificity ranging from iV-glycans to 
chitooligosaccharides. Three proteins clustered under the 
term carbohydrate metabolism were negatively correlated 
with VSH in larval tissue. Ecdysone-inducible gene 
(ImpL3), phosphoglucomutase (Pgm) and gene analogous 
to small peritrophins (Gasp) all process carbohydrates. 
Gasps specifically bind or regulate chitin structure in 
developing embryos, ImpL3 is a lactate dehydrogenase 
induced by the prohormone of 20-hydroxyecdysone that 
regulates insect molting, while Pgm is a glycolytic 
enzyme, a process vital for the biosynthesis of chitin from 
glycogen. Taken together, these observations suggest a 
role for larval chitin biosynthesis and/or structural regu- 
lation in hygienic and VSH behavior. In the proteins 
positively correlated with VSH, gene enrichment identi- 
fied casein kinase II, alpha 1 polypeptide isoform 1 
(Ckllaphal or CK2) and peptidoglycan-recognition pro- 
tein-SA (PGRP-SA). Both are components of the 'cell 
surface receptor linked signaling pathway' ontology and 
effectors of interferon and lipopolysaccharide (LPS) 
macrophage inflammatory signaling [24]. PGRP-SA 
detects Lys-type peptidoglycan (PG) from gram-positive 
bacteria, leading to activation of the Toll receptor path- 
way and, ultimately, to increased expression of antimicro- 
bial peptides [25,26]. The larval inflammatory response 
may serve not only as an individual defense mechanism 
but also as an initiator of social immunity behavior, that 
is, VSH. 

Discussion 

We have described here the discovery of several proteins 
whose expression levels may impact honey bee resistance 
to infestation by the Varroa mite. Natural diversity in 
these behaviors was a prerequisite to this study and we 
observed that the levels of each behavior in any given col- 
ony were not random. As expected, there was a strong 
negative correlation between mite infestation levels and 
HB. At the expression level, several proteins were highly 
significant predictors of HB and mite infestation 
dynamics. Highlighted within these proteins were the 
putative ApoO homolog and a putative Tg. Apolipopro- 
teins are called apolipophorins in insects, and they have 
diverse roles in lipid solubilization and the transport of 
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Table 4 Summary of gene enrichments for antennal proteins correlating with anti-parasitic traits. 


Behavior 


Expression 


Larval Tissue Enriched Gene ontologies 


R24 


t18 


Peptidyl-amino acid modification (caf1, Cg6370), structural constituent of chitin-based cuticle (Cpr65Av, Ccp84Ae), 
membrane enclosed lumen (CG1 140, CG21 18, Adk Caf1) 




116 


Protein modification process (Rep, Art, Pp2A-29B), regulation of neuro transmitter release (Rep, beta-spec), phosphorous 
metabolic process (Pp2A-29B, CG7712, PyK), 


U24 


|18 


Peptidyl-amino acid modification (caf 1 , Cg6370), cell division (Caf1, hts), cell periphery (Hexol, NepYr) 




110 


Protein modification process (Rep, Ugt, eff), regulation of neurotransmitter release (Rep, beta-spec), nervous system 
development (beta-spec, eff, RpL30) 


R48 


|12 


Regulation of embryonic development (Hrb27C, Cg 1 881 1 ), oogenesis (hts, Hrb27C, nudC) 




118 


Muscle tissue development (beta-Spec, alt, wupA), neurological system process (Rep, alt best-spec), cytoskeletal organization 
(Beta-spec, alt, cher, Rpl32, Rpl30, elF3-S10) 


U48 


f13 


Regulation of embryonic development (Hrb27C, CG1881 1), oogenesis (nudC, hts, Hrb27C) 




117 


Muscle tissue development (alt, wupA), ligase activity, forming carbon-nitrogen bonds (pug, CG1315) 


VSH 


-f 18 


Cell surface receptor linked signaling pathway (PGRP-SA, Ckllapla), bristle development (Chllalpha, CG12163) Translation 
elongation (Ef 1 beta, Ef2b) 




415 


carbohydrate metabolic process (lmpL3, Pgm, Gasp) 



R, removed; U, uncapped; VSH, Varroa sensitive hygiene. 



small hydrophobic ligands [27-29]. In innate immunity 
the apolipophorin ApoLp-III stimulates antimicrobial 
activity in the hemolymph, acting as a pattern recogni- 
tion system for LPS and lipoteichoic acid (LTA) [29]. 
Lastly, the strong correlation of Tg with both NDs and 
an increase in the ratio of phoretic mites to brood mites 
suggests that Tg activity could provide a measure of 
resistance to Varroa reproduction. 

V. destructor is an ecto-parasite feeding communally 
and repeatedly on hemolymph of the honey bee through 
a bite wound in the cuticle [30-32]. In insects' innate 
immunity the cuticle provides the first line of defense; 
once breeched, innate defense systems of the haemocoel 
cavity are orchestrated by hemocytes, the fat body and 
hemocoel [33]. Normal wounds heal as hemocytes and 
plasmatocytes exocytose the clotting factors hemolectin 
and Eig71Ee [34]. These molecules and other plasma- 
based factors such as fondue are cross-linked by Tgs in a 
Ca 2+ dependent mechanism to form a primary clot. How- 
ever, V. destructor transmits bio-active compounds that 
prevent healing and allow continued feeding to occur at 
the same wound [35]. In the tick arthropod-mammalian 
ecto-parasitic systems, 18 known bio-active suppressants 
target innate antiseptic defenses, including several 
immune cells types, inflammatory and coagulatory cas- 
cades [36]. In honeybees, the effect V. destructor elicits 
on the immune system is uncertain. Yang and Cox-Foster 
[37] demonstrated that Varroa parasitism increases the 
susceptibility of adult bees to bacterial infection, but no 
major immunosuppressive effects were revealed by tran- 
scriptomic studies on specific immune genes or in global 
analyses [38,39]. More recently a study has reveled that 
salivary secretions from the Varroa mite are able to 
damage hemocyte aggregation in the tomato moth, 
{Lacanobia oleracea) [38] but no known factors of either 



pathogen or host are identified. We report here that ele- 
vated expression of a putative key clotting factor (Tg) is 
found in the larva of Varroa resistant bee colonies. These 
data indicate that honey bees have adapted to Varroa, 
increasing the clotting capacity of hemolymph in order 
to limit mite reproduction. 

While the experiments described here were clearly of 
sufficient power to permit the discovery of some correla- 
tions between protein expression and behavioral traits, 
the variability within such out-bred populations is very 
high. This is likely a significant limitation in fully defining 
the molecular mechanism of something as complex as a 
behavior. Practical limitations in the number of colonies 
that could be sampled and the depth to which the pro- 
teome could be measured across multiple samples were 
inherent problems here, as with any proteomics study. 
Even so, an exploratory approach was seen as an impor- 
tant step in generating new hypotheses in a currently 
poorly understood area of biology. 

It is thought that the speed with which hygienic bees 
respond is driven by a lower limit of olfactory detection of 
the diseased brood odor [40], which is in turn influenced 
by the neuromodulator octopamine [41]. In the antennal 
lobe, octopamine concentration varies between behavioral 
state, being low in nurse bees and high in foragers. 
Juvenile hormone and brood pheromone both modulate 
behavioral responses to octopamine [42] and both are 
involved in several aspects of behavioral maturation, with 
the best-understood system being the transition from 
nurse to forager. This maturation invokes physiological 
changes that are underpinned by increased neural proces- 
sing which is required to interpret complex visual infor- 
mation for flight behavior. Anatomically, expansion of the 
mushroom body neurophil space in the brain and decrease 
in the volume of the olfactory glomeruli of the antennal 
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lobes occurs during this transition [43] . Olfactory sensory 
neurons from the antennae project onto the glomeruli of 
the antennal lobe via the antennal nerve, and olfactory 
information is processed and projected to higher-order 
brain centers such as the mushroom bodies or lateral 
protocerebrum. 

The data presented here indicates that cells (most likely 
neurons with antennal axons) of bees performing rapid 
hygiene express different levels of proteins involved in 
adhesion and vesicle processing (Figure 6a), supporting 
the role of octopamine and maturation as an important 
control of this behavior. The cell adhesion proteins identi- 
fied were all integrin proteins, some of which have been 
reported to regulate synaptic plasticity [44] . Specifically, 
ankerin 2 stabilizes synaptic connections to the spectrin- 



actin cytoskeleton and laminin A, Zasp and Fasl are 
involved in the assembly of functional integrin adhesion 
sites essential for growth cone extension in axon guidance 
during neurogensis [45] . The increased expression of vesi- 
cle sorting proteins in hygienic bees indicates that while 
plasticity may be reduced, antennae of hygienic bees pro- 
vide a strong input into higher brain function. These data 
could be explained by the environment of a hygienic nest 
bee, in which strong brood and queen-based olfactory 
cues are the major sensory inputs for bee development, 
behavior and social cohesion [43] (Figure 6d). Dimorphism 
in neural plasticity has been well characterized in the 
antennae of drones, where the antennal sensory nerves are 
thicker but project into a smaller number of glomeruli 
than in workers [46,47]. This configuration provides 
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Figure 6 Summary model for proteins identified in antennae and larvae from correlation screen. For details see discussion 
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drones with the lower limit of detection for queen phero- 
mone, enabling efficient queen finding during mating 
flights. 

VSH limits mite reproductive success in the brood by 
specifically detecting the presence of a post-ovipositional 
mite. As part of the bees' response, a sensitive adult 
uncaps and re-caps the cell, effectively inhibiting mite 
reproduction [48] . The signal being sensed in this process 
remains unknown, although it peaks between three and 
five days after the cell is initially capped, leading to spec- 
ulation that VSH adult bees respond to temporal fluxes 
in pathology mediated by oviposition, wounding related 
stress responses, infections, and olfactory cues [48]. Cor- 
relation between VSH scores and two proteins encoding 
divergent members of the To/JHBP super-family suggest 
they may be functionally linked to the behavior; To/ 
JHBPs contain a conserved ligand binding domain with 
differing affinities to small lipophilic molecules such as 
JH and the N-terminal signal peptide indicates that they 
are probably secreted into the hemolymph where they 
act as soluble receptors for their ligands [49,50]. In the 
honey bee genome there are eight To/JHBP genes, 
located at two distinct loci, and we see one protein from 
each loci, one positively correlated with VSH and one 
negatively correlated (Figure 6b). Biologically, this separa- 
tion in the genome suggests divergent functions and this 
is further supported by their differential regulation in our 
study. One of these proteins is completely uncharacter- 
ized but in Phormia regina the ortholog of the other To/ 
JHBP is thought to be involved in chemosensation in 
antennal olfaction and taste [49] leading to the attractive 
hypothesis that it is playing a similar role in sensing 
brood. 

That sensory and neuronal processes have a link to dis- 
ease tolerant behavior may be expected but, intriguingly, 
a class of proteins involved in larval cuticle formation/ 
structure also emerged as likely candidates. An arthro- 
pod's cuticle forms the primary physical barrier to the 
environment so while the cuticle plays an obvious role in 
pathogen defense, how it may contribute to social immu- 
nity mechanisms is less clear. Cuticular lipids differ 
between bees depending upon caste and attacks by V. 
destructor can alter the composition in adults and larvae 
[51,52]. The role of the cuticle in social immunity is sup- 
ported by the data presented here, which indicates that 
several proteins involved in forming and maintaining the 
cuticle are significantly correlated with disease tolerance 
behaviors of nurse bees (Figure 6c). 

Conclusions 

Our analysis of tissue proteomes from a large cohort of 
commercial honey bee colonies provides new clues to 
the evolution of biochemical components facilitating 



adaptation to disease. The control of behavior poten- 
tially represents the most complex paradigm in all living 
creatures so its study in natural, outbred systems is 
fraught with many difficulties, explaining the lack of 
coherent mechanisms describing these processes. Honey 
bees live in eusocial colonies and provide a scalable sys- 
tem for the study of developmental social biology and 
the divisions of labor it defines. Our results represent 
indications of molecular mechanisms underlying innate 
and social immunity behaviors in honey bees and build 
upon previous work demonstrating adaption involving 
neural remodeling and odorant recognition. A focused 
investigation of the processes identified here will provide 
an explanation of how host-pathogen interactions drive 
selection to generate disease tolerant colonies. 

Materials and methods 

Reagents 

All chemicals used were of analytical grade or better and 
all solvents were of HPLC-grade or better; all were 
obtained from ThermoFisher-Scientific (St. Waltham, 
MA, USA). 

Honey bee-Varroa populations and physiology 

We established 40 genetically heterogeneous honey bee 
colonies at a research apiary (Grand Forks, BC, Canada) 
in the spring of 2009 by shaking workers into a large 
cage and then portioning them back into single Lang- 
stroth box colonies with nine frames in each. Selected 
queens were then introduced into each new colony with 
initial populations of 1 kg of bees with relatively uni- 
form V. destructor infestation rates, varying among colo- 
nies from 6.2% to 7.6% per 100 adult bees. Colonies 
were allowed to develop for six weeks to allow worker 
populations to turn over and be composed of the intro- 
duced queen's progeny, at which point we evaluated 
each for physiological V. destructor interactions and HB. 
HB was measured as the proportion of sealed brood 
cells uncapped (U) and removed of pupae (R) within 24 
and 48 hours of freeze-killing defined patches brood 
with liquid nitrogen; PH and BI were estimated as 
described [53]. The proportion of uncapped cells 
referred to all cells uncapped by nurse bees including 
those where the pupae had been removed and those 
where the pupae was still present at the time of the 
observation. To estimate the ND V. destructor in each 
colony we counted the number of mites captured on 
screened bottom boards over four 24-hour collection 
periods spanning a period of 10 days. The ND estimates 
were normalized by colony size using the total weight of 
bees to determine the number of bees in each colony. 
VSH was estimated as the production of sexually viable 
female offspring as described [12]. 
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Analysis matrix 

We used triplex dimethylation labeling and generated a 
D-optimal design matrix to group the samples in blocks 
of three and assigned a label to each sample as described 
[20]. A randomized incomplete block design similar to 
what we have used previously was chosen to minimize 
the standard error of the estimate of the colony effect on 
protein expression level [20]. 

Sample collection and protein preparation 

The antennae and larvae from colonies were sampled in 
triplicate. Ten pairs of antennae from nurse bees and 
three fifth instar larvae were removed in situ and frozen 
on dry ice. Larvae were further dissected to remove the 
digestive tracts and free-flowing hemolymph with stabi- 
lity maintained in PBS (50 mM K 2 HP0 4 , 150 mM NaCl, 
pH 7.4) containing complete, EDTA-free protease inhibi- 
tor cocktail (Roche, Mississauga, ON, Canada). Both 
samples were washed three times in PBS and prepared 
using essentially the same method. Both tissues were 
homogenized in 50 mM Tris-HCl, 150 mM NaCl, 1% 
NP-40, 20 mM dithiolthretitol in a Fast Prep bead mill 
with 2.8 mm ceramic beads (MP Biomedicals, Santa Ana, 
CA, USA) using 1 or 3 cycles of 20 s at 6.5 M/s with 
cooling for larvae or antennae, respectively. Tissue lysates 
were clarified at 5,000 relative centrifugal force (rcf) for 
5 min at 4°C before ethanol/sodium acetate precipitation 
[54]. Proteolytic digestion of 20 \i% (antennal) and 50 (ig 
(larval) total protein was carried out as described [54]) 
and samples were labeled by reductive dimethylation 
using formaldehyde isotopologues [55] with slight modi- 
fications [20] . After labeling, each sample was pooled as 
required by the experimental design and each pool was 
separated into five fractions (antennae) using Ci 8 -SCX- 
Ci 8 STAGE tips [56] or into six fractions (larvae) by 
isoelectric focusing with the OFFGEL system (Agilent 
Technologies, Santa Clara, CA, USA) [57]. 

Proteome screen 

Quantitative proteomic datasets were generated for anten- 
nae exactly as described in [20]. For larval tissue, LC-MS 
was done on a 1200 Series nanoflow HPLC system (Agi- 
lent Technologies) interfaced with a chromatin immuno- 
precipitation (CHIP)-cube to a 6520 Q-TOF (Agilent 
Technologies). Peptide separation was performed by 
reversed phase chromatography using a micro-fluidic 
CHIP comprised of an analytical column (75 |im ID, 150 
mm length with a 300 A C 18 stationary phase) and a 160 
nL trap column of the same phase. Peptides were loaded 
in 5% (v/v) acetonitrile, 0.1% (v/v), formic acid at 0.3 \iL/ 
min and then resolved at 0.3 [iL/min for 90 min, during 
which a linear gradient of acetonitrile was created from 
5% to 50% in 0.1% (v/v) formic acid. Mass spectrometry: 



Operating in auto MS/MS acquisition mode, the Q-TOF 
was set up to acquire full scan data over a mass range of 
350 to 2,000 m/z and MS 2 for the six most intense, multi- 
ply-charged ions. Peak lists were created using Spectrum 
Mill extractor specifying fixed modification carboxyamido- 
methylation (C) and triplex dimethyl-mix (K,N-term) 
which accounts for all possible label moieties. Scans were 
merged within +/- 45 s elution time and maximum m/z 
window of 0.5 Da (usually 20 ppm). For database search- 
ing, oxidation (M) was added as a variable modification, 
peptide tolerance was +/- 20 ppm and fragment ion toler- 
ance was 50 ppm, dynamic peak thresholding was 
switched on. Search results from Spectrum Mill were vali- 
dated using autovalidation for protein score >20; (charge 
2, Score >5, % SPI <50, Fwd - Rev Score >1, Rand 1-2 
score >1), (charge 1, Score >6, % SPI <60, Fwd - Rev Score 
>1, Rand 1-2 score >1), (charge 3, Score >5, % SPI <50, 
Fwd - Rev Score >1, Rand 1-2 score >1) and for peptide 
hits (charge 2, Score >12, % SPI <60, Fwd - Rev Score >2, 
Rand 1-2 score >2), (charge 1, Score >8, % SPI <70, Fwd - 
Rev Score >2, Rand 1-2 score >2), (charge 2, Score >12, % 
SPI <60, Fwd - Rev Score >2, Rand 1-2 score >2). Both lar- 
vae and antennae datasets have peptide identification false 
discovery rates (FDR) well below 1% and were compiled 
into protein data arrays as described [20]. All experimental 
design and proteomic results are listed in Additional files 
1, 2, 3 and 4. 

Data availability 

All MS/MS data used in this study have been made 
available in two locations: they have all been deposited 
into the Honey Bee PeptideAtlas [58,59] as processed 
spectra and the raw files themselves are available on our 
FTP site (ftp://foster.chibi.ubc.ca/Downloads/BeeBiomar- 
kers/). 

Statistical analysis 

Logarithms of intensities were normalized by first sub- 
tracting the average of the three measurements in each 
block (for each protein independently) and then center- 
ing and standardizing within each label (across proteins) 
by the median and median absolute deviation [20]. For 
each protein, a linear mixed effects model was used to 
estimate the effect of each predictor variable on the pro- 
tein expression level, adjusting for block and label fac- 
tors. Colony was treated as a random factor to control 
for the three repeated measures within each colony. For 
each predictor variable an estimated effect, standard 
error and P-vaiue was computed for each protein 
response. FDRs (q-values) were computed for the set of 
P-values of a given predictor over all protein response 
variables to adjust for multiple comparisons. All calcula- 
tions were performed in the R statistical language. 
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Gene enrichment analysis 

Analysis was performed using Exploratory Gene Associa- 
tion Networks (EGAN) software [60] with pre collated 
networks for Drosophila melanogaster (Dmel). A. mellifera 
(Amel) gene identifiers were mapped to Dmel orthologs 
using the Round-Up database [61] and any unmapped 
Amel identifiers were assigned functions based on their 
closest homolog in D. melanogaster using BLAST-P, 
resulting in a total of 90% coverage for all antennal pro- 
teins identified. The remaining 10% were dealt with manu- 
ally by drawing information from several sources: 1) honey 
bee genes othologs implicated in immunity [61]; 2) pro- 
teins found significantly regulated in response to bacterial 
infection by Paenibacillus larvae [63]; 3) proteins regu- 
lated in response to V. destructor infestation [18]; and 4) 
proteins specific to colony collapse disorder (CCD) 
affected colonies [64]. Enrichment analysis of proteins 
whose expression levels correlated (P < .05) with behavior 
and their direction of regulation was carried out by 
integrating these nodes with gene ontology nodes for 
component, function and process. For each node, over- 
representation or enrichment analysis was carried out 
employing a standard one-tailed Fisher's exact (hypergeo- 
metric) test using the entire gene list as background. Heat 
maps representing the interaction of important genes and 
their relationships by nodes were generated for both data- 
sets independently. 

Additional material 



generate LC-MS datasets. Each block is shown with colony numbers and 
the output file name generated for each block. Pages 2 to 4 are 
examples of the output from LC-MS and protein identification and 
quantitification result for each peptide used to generate proteins ratios. 
Page 5 is the result of rolling all peptide identification up into unique 
proteins using parsimony rules. 
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Apo: apolipoprotein; Bl: brood infestation; CCD: colony collapse disorder; 
CHIP: chromatin immunoprecipitation; Fasl: fasciclin 1; FDR: false discovery 
rate; Gasp: gene analogous to small peritrophins; HB: hygienic behavior; 
HPLC: high pressure liquid chromatography; JH: juvenile hormone; LC: liquid 
chromatography; LPS: lipopolysaccharide; LTA: lipoteichoic acid; MS: mass 
spectrometry; ND: natural drop; PBS: phosphate buffered saline; Pgm: 
phosphoglucomutase; PG: petidoglycan; PGRP-SA: peptidoglycan recognition 
protein SA; PH: phoretic; PPMC: Pearson's Product movement correlation; 
QTL: quantitative trait loci; RCF: relative centrifugal force; STAGE: stop and go 
extraction; To/JHBP: take-out/juvenile hormone binding protein; VSH: Varroa 
sensitive hygiene; Zasp: Z band alternately spliced PDZ-motif protein. 
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